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ABSTRACT 


The theory of optimum radar detection is well known and is generally implemented in 
expensive ASICs or supercomputers. However, today’s state-of-the-art FPGAs are capable of 
performing relatively complex algorithms and provide the added benefit of being reconfigurable 
with new algorithms or methods on-site. Los Alamos National Laboratory has undertaken the 
goal of developing a receiver that is capable of performing detection and bandwidth estimation of 
pulsed radar systems. It is designed to function in electronic intelligence (ELINT) applications, 
where the goal is to determine the capabilities of threatening systems, such as radars which 
guide aircraft or missiles to targets. 

This thesis addresses methods of pulse detection and bandwidth estimation that are able 
to be implemented on an FPGA. The framework is that which is commonly used in this appli¬ 
cation: a polyphase filter bank subband frequency decomposition of the RF signal, followed by 
statistical detection methods. The optimal fixed-sample-size (FSS) estimator for this subband 
decomposition is shown to be the F-test, based on the output statistics of the filter bank, which 
are found to be chi-squared. An alternative to fixed-sample-size methods, the sequential prob¬ 
ability ratio test (SPRT) is, however, more suited to ELINT due to its ability to adapt the test 
length to the received data. The SPRT is shown to achieve a higher probability of detection 
with approximately 1/5 the required sample size of the FSS method. The complexity of the 
SPRT is equivalent to that of the FSS method, and the statistic that results from the optimal 
SPRT implementation also lends itself to easy calculation of the bandwidth of the signal. 
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CHAPTER 1 


INTRODUCTION 


For years, the theory of optimum radar detection has been known and performed in prac¬ 
tice using analog methods, expensive ASIC implementations, or supercomputers. The recent 
advances in FPGA technology have caused a new paradigm in algorithm development, the re- 
configurable computing architecture. Algorithms coded in the hardware description language 
(HDL) can be loaded onto FPGAs to perform signal processing tasks at a speed that is not pos¬ 
sible today in low-cost floating-point processors. The ability of the FPGA to be reprogrammed 
with a new or updated algorithm allows it a flexibility that is not available in custom-designed 
ASICs. With this reconfigurable flexibility come certain limitations, however. The compu¬ 
tational complexity that is achievable in today’s state-of-the-art FPGAs is a function of the 
number of logic components that can be fit on a chip. Also, the FPGA only functions efficiently 
in a fixed-point architecture, which does not lend itself to complex signal processing algorithms. 
As such, it requires a joint optimization of resources to reduce the computational complexity of 
the algorithms so they will perform well in the environment of reconfigurable computing, while 
maintaining a level of performance that justifies the use of RC technology over more expensive 
ASICs. 

A leader in RC technology, Los Alamos National Laboratory’s Space Engineering Division 
has recently developed an RC platform suitable for signal processing applications and has un¬ 
dertaken the goal of implementing a matched bandwidth radar receiver in FPGA technology. 
The design goal is to implement a computationally efficient radar pulse detection strategy, while 
maintaining as high a level of performance with respect to optimal methods as possible. This 
thesis work provides a piece of the overall design by investigating the best pulse detection strat- 
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egy followed by a bandwidth estimation scheme that can be implemented given the restrictions 
of the RC platform, the RCA-2. 

1.1 Problem Statement 

The particular design problem posed was that commonly incurred in an electronic intel¬ 
ligence (ELINT) application. ELINT is the result of observing signals transmitted by radar 
systems to obtain information about their capabilities [lj. The value of ELINT is that it pro¬ 
vides information about threatening systems, such as radars that guide aircraft or missiles to 
targets. Clearly, ELINT is most useful in situations where some hostility is involved; otherwise 
the information could be obtained directly from the radar user or designer. The design problem 
is difficult in that, for ELINT applications, a priori knowledge of the systems to be detected is 
minimal. Besides the knowledge that the detector should find pulsed radar systems, a minimum 
of other assumptions should be made. The design method can be summarized most effectively 
as a series of steps: 

• Precondition the signal, possibly involving linear or nonlinear transformations. 

• Implement the “best” detector, given the computational complexity limitations of the 
host platform. 

• Upon detection, deduce as much information regarding the detected event as possible (e.g., 
time of arrival (TOA), pulse length, SNR, bandwidth, center frequency, pulse modulation). 

In fact, the optimal solution to this problem is known; however, the goal of this work is to see 
if the optimal solution can be implemented in an FPGA, and whether it is wise to do so. There 
may be assumptions and simplifications that allow near-optimal performance for much reduced 
computational complexity. As always, engineering design is a balance of trade-offs between 
performance and implementation issues. 

1.2 Design Parameters 

Several of the design parameters of the problem are limited by the RC platform itself 
and must be heeded during the design process, outlined in Table 1.1. Other parameters are 
assumptions based on the typical environment of this system. Most of the design parameters 
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Tab 

e 1.1 Design Parameters 


Parameter 

Value 

Unit 

Sampling Rate 

100 

MHz 

Bandwidth 

5-45 

MHz 

Pulse Duration 

.1 - 1000 

/IS 

Pulse SNR 

-20 - 20 

dB 

Center Frequency Agility 
Modulation 

Noise 

1 

NOMOP, LFMOP, PSKMOP, Hopped 
White, Gaussian, Non-stationary 

//s 


are self-explanatory; however, those having to do with the center frequency agility and the 
modulation will be dealt with fully in Chapter 5 and are peculiar to this application. 

1.3 System Overview 

As stated earlier, the design goal is a matched-bandwidth pulse detector. A matched- 
bandwidth pulse detector is one in which pulses are detected in the maximum system bandwidth 
but are then modulated to baseband and filtered with a filter matched to their modulation or 
pulse bandwidth in order to improve the resulting SNR of the pulse. After this process, there 
are many options, including compression, recording, or further processing of the enhanced 
waveform. At the system level, the process is a parallel one. We are allowed to manipulate the 
input RF data in any way, destructive or otherwise, to accomplish detection, as the same RF 
data is input into a digital modulator which translates the pulse center frequency to zero, then 
uses an adaptive decimation scheme to match the bandwidth of the pulse. Figure 1.1 presents 
a system block diagram. The individual components of the system will be discussed more 
fully, beginning with the polyphase filter bank in Chapter 2. Two options for statistical pulse 
detection are presented: the fixed sample size estimators in Chapter 3, and sequential estimators 
in Chapter 4. Following detection, bandwidth and center frequency must be deduced, and are 
the subjects of Chapter 5. 
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4 










CHAPTER 2 


POLYPHASE CHANNELIZER 


The first step in the process of matched bandwidth detection is the subband channelization, 
which is performed to increase the SNR of pulses, as we expect them to have a finite frequency 
support limited to some ratio of the full input bandwidth. Of course, there are many ways that 
one can go about the filtering that produces a subband decomposition, but the theory rests on 
two operations: a prototype filter to achieve the desired subband shape, and the modulation 
method which produces the channels across some portion of the spectrum. In our case, there are 
three main options, a wavelet decomposition, a filter bank based on a discrete cosine transform 
(DCT) modulator, and a filter bank modulation based upon a discrete Fourier transform (DFT). 
For FPGA implementation, we choose to use a DFT-based scheme because there exist many 
cores for commonly used FPGAs to perform both complex and real-valued DFTs. Wavelet and 
DCT methods are possible, and are investigated fully in a paper by Arrowood [2]. 

Given that we are using a DFT as the modulator for our decomposition, a polyphase im¬ 
plementation makes the most sense and provides a multirate architecture that can allow the 
internal clock used for computation to be set at a much lower rate than the input data clock. 
A polyphase filter bank is a simple system, based on two steps; first, the signal is filtered by a 
prototype filter, and then the resulting signal is transformed via the DFT. A third step, which 
is important to the nature of this project, is what to do with the DFT output. The resulting 
complex output has twice as many data points as the input, because each input sample produces 
a real and imaginary part which must be accounted for in the FPGA individually. As such, one 
option is to use the magnitude of the output, thus resulting in a 1:1 data transformation. 

Since the goal of our subband decomposition is detection, we must derive the expected 
statistics of the output of the polyphase filter bank in order to determine the optimal methods 
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Figure 2.1 Polyphase filter bank. 


of detection. To begin, there must be some assumptions made. We will assume the wideband 
RF noise to be Gaussian, white, and nonstationary. The nonstationarity is important for real- 
world applications, but in most derivations, we will assume that over a short time interval a 
random sample is independent and identically distributed (i.i.d.). Our method of derivation is 
simple. Model each process as a transformation, and apply methods of statistical derivation 
appropriate to each type of transform. 

2.1 Reduced Signal Model 

As a first approximation, let us propose a reduced signal model where the pulse is of finite 
time and frequency support, with no modulation either on the pulse envelope or frequency basis. 
Surface search radar employs this type of pulse often; however, it is far from the expected in 
the case of more sophisticated radars, as modulation provides important sensitivity increases. 
Let our received signal be of the form 

r(k) = s(k) + n(k) (2.1) 

where s(k) is the desired pulsed signal and n(k) is the noise, subject to the assumptions above. 
Let us make several assumptions regarding the signal and noise that provide a reduced signal 
model but give us a roadmap for more complicated signal models. Specifically, we assume: 
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• The DFT length is Nf and the polyphase prototype filter is a boxcar filter of length 
N h = N f . 

• The signal pulse is of length N w = Nf. 

• The signal is unmodulated either in pulse envelope or in frequency, but can contain some 
initial phase offset 6 which is modeled as a uniform density on the interval [0, 2n). 

• The base frequency of the signal is centered in one channel of the decomposition. 

Given the above assumptions, s(k) is given by 

s(k) = Acos(uj(k) + J) (2.2) 

where A is the amplitude, S is the initial phase, and w(k) is of the form 

,2*Jfc£ 

u(k) — e N i (2.3) 

In the absence of noise, the received signal is r(k) = s(k). We model the DFT as a linear 
transformation given by 


Wo,i 


W 0<Nf 


where Wi /. is 


W N . o 


W N , N} 


i 2?r kl 

W iJe = e N i 


If we apply a sequence of linear transformations, W for the DFT and H for the polyphase 
filtering, the real and imaginary parts of the transformed signal are 

X s = Re{ WHs) (2.5) 


Y s = 7 m (WHs) 
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where s = [s(0), s(l),. .. ,s(Nf)] is the signal in vector form. Since we have already assumed 
the prototype filter to be boxcar, H is the Nf x Nf identity matrix. By trigonometric identities, 
we represent s(k) as 

s <*> = I 

Applying the transformation, we recognize that for f = l, the inner product of W and s will 
be a constant, but for f ^ l, it will be zero by orthogonality. Thus, we have 

X s = Re(Ae^ 5 ) = Acos(5) (2.8) 

and 

Y s = Im(Ae' s ) = Asin(<5) (2.9) 

for the real and imaginary outputs. The A absorbs the constant due to the DFT gain and 
polyphase filtering into the original signal amplitude A. Stochastically, X s (l) and Y s (l) are 
characterized by degenerate distributions, as they place all probability at a single point. 

We have derived the statistics of a received pulse in the absence of noise. Now, we look 
at the results of a received signal of noise alone, r(k) = n(k). The signal n(k) is a random 
i.i.d. sample of size Nf, with each marginal distribution Gaussian. The joint pdf of the random 
sample is multivariate Gaussian of dimension Nf. A theorem from multivariate distribution 
theory will be subsequently helpful. 

Theorem 2.1 (Linear Transformations of Multivariate Normals) Let X = [X x , X 2 ,..., X n ] 

be a random sample with a joint pdf that is multivariate normal and is distributed as 

X ~ N n (p,Y,) ( 2 . 10 ) 

where the mean vector is p and the covariance matrix is 2. A random variable Y obtained by 
a linear transformation of the random sample X, Y = LX is distributed as 

Y~jV n (L/z,LSL') (2.H) 

where ' denotes transposition. 

Because the DFT is a linear transformation, we can use the theorem to find the distribution of 
the real and imaginary outputs due to noise: 

X n = Ue(WHn) ( 2 . 12 ) 


• 2jvk_l 

e 3 N f 


• 2rr kf 

e jS + e 3 N f 


e* s 


(2.7) 
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and 


Y n = Jm(WHn) (2.13) 

where W and H are as before. Thus, the distributions of X and Y are given by 

X n ~ N(Re( W)/x, f?e(W)Sf?e(W)') (2.14) 

and 

Y n ~ N(Im(W)n, fm(W)S/m(W)') (2.15) 

Since the antenna is ac coupled to the ADC, the noise is zero mean, and the covariance matrix 
is a 2 1 by the i.i.d. assumption. Furthermore, since W is an orthogonal matrix, 

X n ~ N{ 0, a 2 Re{ WW')) = N( 0, a 2 ) (2.16) 

and 

Y n - N{ 0, a 2 /m(WW')) = N{0, a 2 ). (2.17) 

To find the distribution of the received signal of the original reduced model form r(k ) = 

s{k) + n{k), we only need to apply the individual results from the signal and the noise alone. 

Theorem 2.2 (Convolution Formula) Let Xi and X 2 be independent random variables with 
probability density functions f x (x) and f 2 (x), respectively. Let there be a transformation, Y x = 
X\ + X 2 and Y 2 = X 2 . Applying a change of variable technique, the joint pdf ofY\ and Y 2 is 

9 {yuV 2 ) = fi{yi -y 2 )/ 2 (y 2 ) (2.18) 

and the marginal pdf of Y\ = X\ + X 2 is given by 

/ '00 

fi(yi - V 2 )f 2 (y 2 ) dy 2 = fi(xi) * f 2 {x 2 ). (2.19) 

-OO 

Thus, the pdf of the sum of independent random variables is the convolution of the constituent 
pdf’s. 

By Theorem 2.2, X and Y are distributed as the convolution of the distributions of s(k) and 
n{k), 

X s + n ~ -j=r- e * S(Acos(S)) = N(Acos( 5),a 2 ) (2.20) 
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and 


Ys+n ~ V^rro- 6 * <*(-Asin(<J)) = iV(ism(<5),cr 2 ) (2.21) 

so X and Y are distributed as normal random variables with equal variance with mean depen¬ 
dent on the initial phase of the signal. 

To complete the derivation of the statistics for the reduced model, we must now apply the 
final step of magnitude conversion. In general, the magnitude of the DFT would be calculated 
as Z = VX 2 + Y 2 , but since square roots are very computationally complex in the FPGA, we 
choose to use the squared magnitude Z = X 2 + Y 2 instead. From distribution theory, we know 
that for squared normal variables, we expect a chi-squared distribution, and in fact, we expect 
two distributions in our case. For noise alone, 


and for the signal plus noise case, 


Z ~ a 2 xl 


( 2 . 22 ) 


21 


Jlj2 


X2 (A) 


(2.23) 


where Xn is the central chi-squared distribution with n degrees of freedom, and x 2 {\) is the 
noncentral chi-squared distribution with n degrees of freedom and noncentrality parameter A. 
In our case, the noncentrality parameter is 


A = 


Acos(<5) 


+ 


Asin(<i)\ A 2 


(2.24) 


Thus, we expect chi-squared random variables at the output of our polyphase filter bank. 


2.2 General Signal Model 

Admittedly, the reduced signal model is not of much use for signals of interest in ELINT 
applications, but it gives us a basis for extrapolating to a more general signal model that encom¬ 
passes such variables as pulse modulation, phase modulation, and pulse-length considerations. 
Of the assumptions of Section 2.1, we retain our assumptions about the noise only, as we will 
throughout this entire thesis. However, we now allow both phase and pulse modulation of our 
signal 

s(k) = A(k) cos(u >(k) + d(k) + 6) (2.25) 
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where A(k) is the pulse modulation, and 6{k) is the phase modulation. 

It is excessive to repeat the above derivation of the statistics as we can use the results to 
make general inferences. The most simple case is that of only pulse envelope modulation, as that 
would imply that the noncentrality parameter is no longer constant but would be a function of 
time. In the case of phase modulation, the inner product of the DFT matrix W with the signal 
would no longer be constant in one channel, but would be constant in any number of channels, 
and also a function of time. Both of these results do not change the expected statistics; we 
retain a chi-squared distribution, but it becomes nonstationary. 


2.3 Prototype Filter and Pulse Lengths 

In the reduced signal model, we assumed a boxcar filter of length Nr, the DFT length. 
However, this resulted in no prototype filter at all. In general, we would like a prototype filter 
that achieves the desired out-of-band rejection characteristics. The filter prototype will usually 
need some overlap of the DFT to accomplish this, and we assumed above that the pulse length 
was equal to the filter length. For pulses longer than the filter length, this is not a problem, 
but for any pulse shorter than the filter length we expect to see a reduction in SNR as the 
filtering operation spreads the concentrated energy of the pulse to the full length of the filter. 
In an optimal situation, since the pulse length is unknown, we would run multiple filter banks 
in parallel, each with a different prototype filter sized to optimally filter a pulse of expected 
length. However, FPGA implementation prohibits this excessive computation. Thus, we must 
choose with care the filter length in order to preserve SNR for short pulses, while providing 
sufficient channel subband characteristics after the transformations. 

An important function of the prototype filter is to shape each subband channel, but equally 
important is the effect of the filtering on the statistics between and within the channels. Our 
detection methods will rely on statistical inference, or hypothesis testing. In many methods, the 
statistical independence of the samples is important in determining multivariate distributions. 
For our purposes, we would like to limit the correlation of the samples as much as possible, in 
order to declare them independent to the tolerance of our tests. We can use the methods of 
Section 2.1 and particularly the theorem of transformations of multivariate normals to find the 
expected correlation of the output variables Z. For the case of a nonboxcar prototype filter of 
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any length, the distribution of the real and imaginary components in each channel is 


X ~ Af((WH)ylcos(£),a 2 (WH)(WH)') (2.26) 

and 

Y ~ iV((WH)j4sin(<5),a 2 (WH)(WH)') (2.27) 

If we work with the covariance matrix £ = <t 2 (WHH'W'), we can find the correlation coeffi¬ 
cients both between and within each channel for any prototype filter and DFT combination. 


2.3.1 Intrachannel correlation 


By intrachannel correlation, we mean the correlation of the samples between the differ¬ 
ent channels of the polyphase output, before the magnitude conversion. For example, how is 
sample k of the eighth channel correlated with that of its neighbor in the seventh or ninth 
channels? This becomes important in bandwidth estimation later, where we would like sta¬ 
tistical independence in order to infer whether there is significant energy in any channel. For 
the covariance matrix S above, statistical independence implies orthogonal transformations. 
Because the DFT matrix W is necessarily orthogonal, if we limit the polyphase transformation 
of the prototype filter H to be orthogonal we require strict statistical independence between 
channels. Specifically, we guarantee that the corellation coefficient given by 


_ Cov(Zt, Zj) 

(JiCTj 

is zero in every other channel for aligned data. 


(2.28) 


2.3.2 Interchannel correlation 

Similarly to intrachannel correlation, by mferchannel correlation we mean the correlation 
between different samples in the same channel. For example, how are samples k and k - 1 of 
the eighth channel correlated? This becomes important in the pulse detection stage, as we will 
employ a channelized detector for each channel, independent of the other channels. For this 
analysis, we must generalize H so that the linear transformation gives us more than one output 
sample in each channel for comparison. Here, we will become more specific to this project for 
ease of notation in terms of specifying a prototype filter and DFT length, but the results can 
be readily generalized to any length by following a similar methodology. For this project, the 
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DFT length is 32, and a suitable prototype filter length is 128. Let us define the prototype 
filter polyphase transformation matrix to be H of size 32 x 128, and partition that matrix into 
four submatrices each of size 32 x 32 as in 

H = [ Hi H 2 H 3 H 4 ] 

Now, we desire to yield more than one output sample in each channel. If we construct a new 
matrix H* as in 

Hi H 2 H 3 H 4 0 • • • 0 0 0 

0 Hi H 2 H 3 H 4 0 0 0 

0 0 Hi H 2 H 3 H 4 0 0 

and a new DFT matrix W* as in 

WO 0 0-0 

0 W 0 0 - 0 

w* = 

0 0 W 0 --- 0 

where 0 denotes a matrix of zeros of size 32 x 32, then we can form a new linear transformation 
that yields multiple output samples in each channel. This new transformation is given by 

X = Re{ W*H*r) (2.29) 

Y = 7m(W*H*r) (2.30) 

where r is the received signal in the presence of noise. We can now compute the covariance 
matrix E* = a 2 (W*H*H*'W*') and use it to find the intrachannel correlation coefficients. 

2.3.3 Correlation results 

Applying the methods above, we can find the correlation due to any prototype filter. We 
begin with the filter designed specifically for this project by Joseph Arrowood, which is a 
length-128 non-orthogonal FIR filter [2]. By orthogonalizing this filter, we can compare the 
results of the orthogonal and non-orthogonal cases to determine the costs and benefits of each 
method. In Figure 2.2, we can see that the nonorthogonal version of the prototype filter has 
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Nor>-OrthogooaJ Prototype Fitter 


Orthogonal Prototype Fitter 




(a) Nonorthogonal Filter (b) Orthogonal Filter 

Figure 2.2 Normalized frequency response of a non-orthogonal filter designed by J. Arrowood, 
and an orthogonalized version. 

much lower sidebands, resulting in a high out-of-band energy rejection in adjacent channels. 
This is important in bandwidth measurement applications, as we would like to keep a strict 
channelization of the energy without excess leakage. 

Using these two prototype filters, we can now compare the interchannel and intrachannel 
correlation coefficients. From the covariance matrix derived in the above section, we can com¬ 
pute the correlation coefficient p with Equation (2.28). The correlation coefficient is bounded 
by —1 > p > 1. Statistical independence implies p — 0, and the closer p is to zero, the more 
“independent” the respective samples are. We put “independent” in quotation marks because 
for p ^ 0 the samples are necessarily dependent, but we may treat them as approximately 
independent for p near zero. In Figure 2.3(a),(b), the interchannel correlation coefficient is 
necessarily zero, as we would expect, in the case of the orthogonal prototype filter, but is 
about .2 for the non-orthogonal filter. While this is higher than we would like, it remains low 
enough that we may treat them as approximately independent as a trade-off for high sideband 
attenuation, which is more important than statistical independence between the channels. In 
Figure 2.3(c),(d), the intrachannel correlation coefficient is similar for both the non-orthogonal 
and orthogonal prototypes. This is to be expected due to the nature of the polyphase filtering. 
For any overlap where the filter length is greater than the DFT length, we have a projection 
operator that projects the original vector onto a smaller subspace. This data reduction im¬ 
plies that some combination of the data is taking place, resulting in correlation of the output 
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samples. However, we notice that in both cases, the correlation coefficient rapidly approaches 
zero away from the current sample, and is still small for the adjacent sample. Thus, we can 
safely approximate statistical independence, especially when p « .03 as in this case for the non- 
orthogonal filter. Prom this analysis, the obvious choice is to proceed with a non-orthogonal 
prototype filter due to its small intrachannel correlation and high sideband attenuation. 
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-9 -« -7 -6 -S -4 -3 -2 -1 0 -9 -8 -7 -6 -5 -4 -3 -2 -1 0 


Relative Sample Number Relative Sample Number 

(c) Non-Orthogonal (d) Orthogonal 

Figure 2.3 (a),(b) Interchannel correlation coefficients in each of 15 channels for both the 
imaginary and real parts of each channel. (c),(d) Intrachannel correlation coefficients based on 
relative sample number, where 0 is the current sample number and previous sample n um bers 
appear as -k for the real and imaginary parts of a representative channel. 
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CHAPTER 3 


FIXED-SAMPLE-SIZE DETECTION METHODS 


Traditionally, this type of subband detection of radar pulses has been accomplished with 
fixed-sample-size (FSS) statistical hypothesis tests. The optimal test is determined by the 
statistics after the linear transformation that is chosen for the subbanding. Thus, for our case, 
the optimal test for the output of a polyphase filter bank is necessarily much different from the 
optimal test for the output of a cosine modulated filter bank, or a wavelet decomposition. As 
such, we refer to the “optimal” solution as it applies to this particular subband decomposition. 
In general, there is a truly optimal solution to the radar detection problem utilizing matched 
filters and multiple independent filter banks for each expected pulse length. However, this is 
not implementable in today’s FPGAs as it is too computationally complex. We will derive 
and characterize the optimal fixed-sample-size test for the particular case of this subband de¬ 
composition. In fact, for this project, the goal of this work was to determine what exactly is 
the optimal solution given this subband decomposition, which is a standard practice in ELINT 
applications. 


3.1 Optimal Fixed-Sample-Size Detector 

Given the statistics present at the output of the polyphase filter bank derived in Chapter 2, 
we can derive the optimal test for the case of signal plus noise over purely noise. For these 
detectors, we will assume that the detector works in a single subband, independently of all other 
channels. For a 32-band decomposition, we would generate 15 usable channels, neglecting the 
dc and Nyquist bands that are subject to fixed-precision filtering anomalies [2]. Thus, in our 
detection scheme, there would be 15 FSS detectors running independently in each channel. In 
general, more complicated methods of alarming on multiple channels at once are possible, and 
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these will be discussed during bandwidth estimation, as they are more effective and appropriate 
in that setting. 

Prom Chapter 2, we know that the statistics present at the output of the channelizer are 
given by 


Z ~ a 2 xl (Noise Only Case) (3.1) 

and 

Z~cr 2 x 2 2 (A) (Signal Plus Noise Case) (3.2) 

where we define Z as the output of the polyphase filter bank subject to magnitude conversion. 
As we saw before, we will treat these Z as statistically independent subject to the considerations 
of the prototype filter maintaining a small correlation coefficient. 

Our first problem in finding an optimal solution is to find a statistic that represents the 
data well, and separates it into two classes based on its properties. The defining characteristic 
of our two cases of signal and signal plus noise is the noncentrality parameter A. This is directly 
equivalent to a shift in mean, as will be shown later when we derive a consistent estimator of the 
noncentrality parameter for continuous interference. Thus, our statistic should be a sufficient 
estimator for A. Furthermore, we require our tests to be noise-riding. Since we cannot assume 
stationarity, we must directly calculate the parameters of the chi-squared distribution from 
the data, in general using maximum likelihood (ML) estimators. We now present formally a 
method for achieving maximal test power over the class of estimators of a chosen size. 

Theorem 3.1 (Neyman-Pearson Theorem) Let X x , X 2 ,..., X n , where n is a fixed positive 
integer, denote a random sample from a distribution that has PDF f(x ; 9). Then the joint PDF 
of X\,X 2 , ■■■,X n is 


L{0\xi,X2,...,x n ) = f(xi;6)f(x2;9)---f(x n ;9). (3.3) 

Let 6' and 6" be distinct fixed values of 9 to that = [9:0 = 9', 9"], and let k be a positive 
number. Let C be a subset of the sample space such that 

( a ) £’Z’,x n J) ^ k f° r each p° int ( x u x 2 ,...,x n )ec 

(b) L(r ; ;xiS’..’’xn) ^ k f° r each P° int ( X U X 2, ...,X n )eC* 
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Channel 



Time 


Figure 3.1 Graphical interpretation of the choice of random samples for testing hypothesis H 0 
against H\. 

(c) a = Pr[(X \, X 2 ,, X n ) G G\ Hq] 

Then C is a best critical region of size a for testing the simple hypothesis Hq : 6 = 9' against 
the alternative simple hypothesis H\ : 0 = 9" [3]. 

The Neyman-Pearson theorem states that we can form our test of a simple hypothesis on one 
parameter of a given distribution if we can construct a likelihood ratio L(9\x\,X 2 ,... ,x n ), 
based on a sufficient statistic of a set of i.i.d. random variables. Furthermore, if we can find 
a constant k such that the size of the test is a, and the constant divides the hypothesis space 
into mutually exclusive regions, C and C*, then the test is a best test of the hypothesis if the 
size of C* is as great as with any other test. Intuitively, we may think of the size of the critical 
region C as the false alarm probability a , and the size of the complementary critical region C* 
as the probability of detection (3, or the power of the test. Thus, Neyman-Pearson gives us a 
way of maximizing the detection probability for a given false alarm rate (FAR). 

Suppose now that we base our statistic on a sum of i.i.d. random variables. Ideally, we 
can compare two random samples of the channelized data and look for a significant difference 
between their noncentrality parameters. As in Figure 3.1, it makes sense to take independent 
observations of data; one will form the basis of the likelihood ratio for the noise-only hypothesis, 
and a second random sampling that will be tested against the first observation. We define 
mathematically two observation vectors: Zn 0 = [Z 0 , Z 1} ..., Z No ] of size N 0 and a second 
observation Zp^ = [Zn 0+t , Zjv 0 +r+i> • .., Zn 0+7 . + ^ 1 ] of size N x where r is a constant chosen to 
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separate the observations by the average rise time allowed by the input bandwidth. By the 
Neyman-Pearson theorem, we construct the likelihood ratio as the joint pdf of the observations 
Z N 1 over Zn 0 which are each i.i.d. distributed chi-square variables. Thus, the null hypothesis 
H 0 is that the noncentrality parameter Ao estimated from Zn 0 is equal to the noncentrality 
parameter Ai estimated from the observations Zjs^. Thus, the statistic 


represents the likelihood ratio 


No+r+Ni 

77T ^ 

i—No-\-r 

~~Vo 

M S Zj 


(3.4) 


Z<(Ai, ^jVo+r; ^ATp-f-r-fl, . . • , ^jVo+r+jVt) 

L(\ 0 ;Zo,Zi,...,Zn 0 ) 


It is worth noting that this method is not optimal under all circumstances. In Chapter 2, we 
found the correlation coefficient to be small due to the filtering operations, but this is under the 
assumption that the input is i.i.d. Under the H 0 hypothesis, the observations Zn 0 and Zn x are 
independent, and the likelihood ratio constructed from the sufficient statistic of the sum of the 
observed data is optimal. However, under the Hi hypothesis, there may be dependence due to 
the received signal because it is not necessarily a process which is random or independent. In 
these circumstances, improved performance can be obtained by coherently summing the energy 
and then squaring the result. However, in cases such as the general signal model of Chapter 2 
in which we allow a center frequency of the received pulse to be off the center of any band of the 
decomposition and we allow phase modulation, or possibly phase drift due to Doppler effects, 
this sufficient statistic is near optimal since we can no longer coherently sum the energy across 
DFT blocks, and must resort to incoherent averaging to achieve the best performance. 


3.2 Detection without Interference 

We will now show how the general results of the likelihood ratio above can be categorized 
into special cases that yield significant results. We have formed our likelihood ratio from the 
Neyman-Pearson theorem, but to achieve a best test we must find some way to calculate the 
constant k such that we produce a critical region of size a, while maximizing the power of our 




Consider the special case of Ao = 0. This could be considered our standard operating mode 
under the assumptions made in Chapter 2 as we assumed that the noise was zero mean. Thus, 
under the null hypothesis, we are testing Ai = Ao = 0, and wish to reject this hypothesis only 
at a significance level a. Under these assumptions with the null hypothesis valid, our likelihood 
ratio is distributed as 


-rr X2N1 /Nl c 


(3.6) 


where F„ ltl , 2 is the F distribution with v\ and 1/2 degrees of freedom [4]. The F pdf is given by 


If{x) = 


r(^) (v A " 172 


r(^)r(f) 


"2 ) 


(! + S*) 


2 


for x > 0 


(3.7) 


Thus, since the distribution of the likelihood ratio is known, we can use it to calculate the 
constant 


k = F2N u 2N 2 -,a (3.8) 

where F2Ni,2N 2 -,a is calculated from the inverse F cdf such that condition (c) of the Neyman- 
Pearson theorem is valid. These F values are easily found in standard tabulations and are also 
a part of many mathematics packages including MATLAB, S-plus, SAS, and in C or Fortran 
routines. A test of this form is commonly called an F-test and forms the basis of much statistical 
theory. 


3.3 Optimal Properties of the F-test 

This section is primarily interested in explaining in a mathematical sense why the F-test 
has desirable properties for FSS detection given the statistics we expect at the output of the 
polyphase filter bank. Let us denote the hypothesis Ho such that Ai = 0 under the assumption 
that there is no interference and Ao = 0. We denote the n-dimensional sample space V n as 
the vector of observations Z = [Zi, Z 2 , ■ ■ ■, Z n \. The choice of a test of H 0 is equivalent to the 
choice of a Borel set C in V n , rejecting Ho if and only if the observed point Z falls into C*.[5] 
The power of the test is given by 

/3(A U C*)= [ Pa(Z) dZ (3.9) 

which is a function of the noncentrality parameter and the critical region C* only. 
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Definition 3.1 (Uniformly Most Powerful Region) C* is a uniformly most powerful (UMP) 
critical region of a given class A if C* € A, and if, for any C*’ e A, P{X\,C*) > [3]. 

We know that the power of the test depends on Ai only through the intermediary of the means 
of the Gaussian variables X and Y of Chapter 2 which are calculated from the amplitude of 
the pulse. The multivariate noncentrality parameter is = o' 2 A' 2 , where A' is any 

nonnegative constant. This is equivalent to saying that the power is constant on surfaces 

n 

Bv= Ai = <r 2 A' 2 (3.10) 

k —1 

The surface B\> is a hypercylinder whose base is a spherical cone of one nappe [5]. 

Theorem 3.2 (Uniformly Most Powerful F-test) Among all critical regions C of size a 
with the property that the power depends only on Ai through the intermediary of a 2 A' 2 , C is 
UMP [5]. 

We now introduce a space of n-dimensional spheres 

n 

S(A' u A' 2 , ..., A'n , A') = A: | £ A\ = a 2 A' 2 ; A\ = A u A' z = A 2 ,..., A! n = A n (3.11) 

fc=i 

where the constants A! i,... , A' n , A' are chosen with only the restriction that A' > 0 and in our 
case A'i > 0. We define the average power over the sphere S under the critical region C to be 
the integral of the power /3(Ai, C*) over S. The cylinder B can be expressed as a union of these 
spheres [5]. 

Theorem 3.3 (Average Power Maximization) In the class of all similar regions of fixed 
size a for testing Ho, C* maximizes the average power on every sphere S{A'\,A' 2 ,..., A' n , A') 
[5]. 

The F-test is thus optimal at maximizing the average power of any test for all alternatives 
Ai > Ao, resulting in a uniformly most powerful test when we are interested “uniformly” in all 
alternatives, represented by a uniform weighting of the spheres. In a case where we are interested 
possibly in one alternative over other possible alternatives, the F-test may be suboptimal and 
other tests may be superior. 

Thus, we find that the F-test has desirable properties for our application. It maximizes 
the average probability of detection for all SNRs, when we are interested uniformly in all 
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Non-centrality Parameter tor a NOMOP Pulse 



Figure 3.2 Noncentrality parameter Ai as a function of wideband SNR for a NOMOP pulse. 

alternatives, implying that we are interested in detecting a 3-dB pulse as much as say, a 20-dB 
pulse. In the event that there is a nonuniform preference for detection, the F-test may be 
suboptimal. 

3.4 Performance of the F-test 

The performance of the F-test can be characterized by the power with which it rejects the 
hypothesis of noise only, H 0 , or equivalently by its probability of detection over the range of 
pulse SNRs. To calculate the power function, we must find the distribution of the statistic T 
under the alternative hypothesis H\ signal present. In this case, we know from our derivation 
of Chapter 2 that in the presence of signal, the joint pdf of the observation vector Z Nl will 
be noncentral chi-squared with some noncentrality parameter \\ > 0. In general, we found 
this noncentrality parameter to be a function of the prototype filter H, so we must calculate it 
based upon its particular gain and other characteristics. Using the non-orthogonal prototype 
filter discussed in Chapter 2, we simulate the range of SNRs relevant for this project, -20-dB 
to 20-dB, and find the mean noncentrality parameter for a NOMOP pulse. In Figure 3.2, we 
see the noncentrality parameter as a function of SNR, which ranges from .07 for a -20-dB pulse 
to 60.5 for a 20-dB pulse. In this case, we have assumed a NOMOP pulse, but the results 
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generalize easily to a case where the energy is distributed among more than one channel due to 
modulation. In these cases, the apparent SNR in each channel may be less than the wideband 
SNR, but the noncentrality parameter will reflect the apparent channel SNR, and the power 
associated with the F-test at that SNR will apply similarly. 

To find the power of the test, we note that the distribution of the test statistic T is 


xInq/No 




(3.12) 


where A) is the noncentral F distribution with u\ and v 2 degrees of freedom with noncen¬ 

trality parameter A. The noncentrality parameter as well as the degrees of freedom are heavily 
influenced by the lengths of the observation intervals for Zn 0 and Zn, due to a property of 
addition of noncentral chi-squared variables. 


Property 3.1 (Addition of x 2 Random Variables) Let X x and X 2 be independent ran¬ 
dom variables distributed as Xm(^l) an< ^ Xn(^), respectively. Then the sum Y = X\ + X 2 is 
distributed as [4] 


Y~xl +n (5 l + 5 2 ) 


(3.13) 


For the hypothesis H\, signal present, this is the distribution of the statistic T. Thus, the 
power of the test, or the probability of rejecting the hypothesis of signal absent when it is false, 
is equal to the integral 

roG 

$(Ai ,C*)= / pf'{x) dx (3.14) 

Jk 

which is the area under the noncentral F pdf over the critical region C*. This critical region is 
the complement of the space C which is the region for accepting Hq when it is true. In general, 
the power of the test is dependent not only on the noncentrality parameter but also the critical 
region specified by the size of the test a, as well as the lengths of the observation samples N 0 
and N\. In Figure 3.3 we calculate the power of the test for different sizes ranging from an a 
of 10 -2 to 10~ 5 and observation intervals N\ = 4 and Nq = 64. 


3.5 Observation Interval Lengths 

From the noncentral F distribution under the hypothesis Hi we see that both the degrees 
of freedom as well as the noncentrality parameter influence the distribution of the statistic. 
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Probability of Detection for Various a with the F Test 



Figure 3.3 Probability of detection (power) of the F-test for various test sizes a and observation 
lengths N\ = 4 and Nq = 64. 


Intuitively, we expect that longer window lengths improve detection probability as they lower 
the variance in the statistic F. In Figure 3.4, we see that as the window length N\ gets longer, 
the probability of detection improves at low SNR. 

However, care must be taken not to overestimate the observation interval length JVi, as this 
will result in a reduction in the apparent SNR of a pulse shorter than N \. In such cases, the 
apparent noncentrality parameter A will be 



(3.15) 


where L is the true length of the pulse. An example is given in Figure 3.5, where a 10-dB, 
1.28 fi s pulse is overestimated by up to 10 times. Notice that even overestimating the window 
length by a factor of 2 appears as less than half the true SNR to the detector. In all cases, N\ 
must be shorter than the observation interval Zjsj 0 in order to keep the number of degrees of 
freedom sufficient in the denominator of the F-test to achieve reliable results. 
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Pulse SNR (dB) 

Figure 3.4 Probability of detection (power) for various length of the Zn x observation interval 
with No = 64. 


3.6 Detection in the Presence of Interference 


In many practical instances the detector may have to perform in the presence of interference 
such as that from continuous wave (CW) carriers due to civilian communications or other 
radars that are not interesting for ELINT purposes and are preferably ignored. In general, this 
interference takes the form of the pulses we would like to detect, such as in the general signal 
model of Chapter 2, but are of longer duration. In this case, to calculate the constant k and 
the power of the associated test, we must again find the distribution of the statistic T. In this 
case, T is the ratio of two noncentral chi-squared variables and is distributed as 


T 


X 2 2 No (N 0 \ 0 )/No 


= R 


2Ni ,2No 


(Ai.Ao) 


(3.16) 


where F” l U2 ( Ai, A 2 ) is the doubly noncentral F distribution with degrees of freedom v\ and z/ 2 
with noncentrality parameters Ai and A 2 . Essentially, the F tests previously were generalizations 
of this distribution with one or both noncentrality parameters equal to zero. The doubly 
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Apparent SNR Due to Overestimated Observation Interval Length 



1.28 2.56 3.84 5.12 6.4 7.68 8.96 1 0.24 11.52 1 2.8 

Observation Interval Length (n sec) 

Figure 3.5 Apparent SNR of a 10-dB pulse 1.28 /zs long caused by overestimating the obser¬ 
vation interval length up to 10 times. 


noncentral F PDF is given by 


Pf»{x) — e 


where 


— i(Ai+A2) 


'E'£E(\ 1 ,r)E(\ . 2 , t 


r=0 5=0 


r(^+r + s ) 

r& + r)r(? + 2) U 


E(\,m) = 


r(m + 1 ) 


tl+^2. +r+s 

( 1 +g*) 

(3.17) 


(3.18) 


valid for x > 0 and Ai > Ao [6],[7]. In this case, under both hypotheses, the statistic is 
distributed as F", with Ai = Ao under Hq, and Ai > Ao under H\. To find the constant k 
in the Neyman-Pearson theorem, we note that a property of F" is that for Ai = Ao, as under 
the null hypothesis, F" reduces to the noncentral F distribution F\ allowing easy calculation 
of k via mathematics software, as in the case of the F-test. However, for the general case 
of an unknown SNR interferer, the noncentrality parameter Ao must be calculated from the 
observations of Zn 0 . Given the moment generating function of x„ 2 (A), A) it can be 

shown that the sample mean is 


dM(t ; i/, A) 


= v + A 


(3.19) 
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(a) Mean Noncentrality Parameter (b) Power 

Figure 3.6 Mean noncentrality parameter and probability of detection (power) of the noncen¬ 
tral F-test for various test sizes a and observation lengths JVi = 4 and No = 64 in the presence 
of a 2-dB CW interferer. 


giving a consistent estimator for Ao as Ao = X - 2. This test is commonly referred to as the 
noncentral F-test. As in the case of no interference, the power of the test is calculated as 

roc 

/3(Ai ,C*)= / p F n(x)dx (3.20) 

Jk 

using numerical integration. In practice, the numerical simulation of Equation (3.17) is sensitive 
to both the magnitude of Ai as well as to the upper limit of summation over r and s. For large 
Ai it will not converge to a meaningful distribution, and an appropriate upper bound for r and 
s is approximately 30 to 50. 

The power calculation is complicated through the possible constructive or destructive inter¬ 
ference of the interferer with the desired pulse. Under the assumption that the initial phases 
of both the interferer and the pulse are independently uniform on the interval [0, 2ir), it is 
possible to achieve anywhere from totally constructive to totally destructive interference. As 
such, unlike the previous cases, even a NOMOP pulse will have a noncentrality parameter Xi(t) 
that is a function of time. As an example, simulation was used to find the mean noncentrality 
parameter of a 3 ps pulse in the presence of a 2-dB interferer when both the interferer and 
desired pulse have unknown initial phase. In Figure 3.6, the power of the noncentral F test is 
calculated using simulation. The values are only approximate and are subject to the variance of 
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numerical integration and the fluctuation of the estimated mean of Ai. However, it is instruc¬ 
tive in that it shows how detection is severely influenced by interferers, as desensitization of the 
receiver results. Even though the noncentrality parameter for a 20-dB desired pulse is large, 
the noncentral F-test performs poorly because of the nonlinear relationship of the noncentrality 
parameter to the shape of the distribution. Even at large SNR, we cannot achieve sufficient 
separation between the F" PDF under the H\ hypothesis and the F' PDF of the Hq hypothesis. 

3.7 FPGA Implementation 

For the optimal FSS tests presented in this chapter, FPGA implementation is relatively 
efficient. In the case of the F-test, a single value of k can be stored and is applicable to all tests, 
but in the case of an interferer, there will be a critical value k for every SNR possible for the 
interferer. This is not a problem for FPGA implementation as the noncentrality parameter is 
calculated from the sample mean and the corresponding critical value is fetched from a lookup 
table. For the FSS detectors presented here, the operation count yields Nq + N\ additions, two 
multiplications, and two divisions. The divisions are most easily accomplished via bit shifting 
when Nq and N\ are of the form 2 n . 
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CHAPTER 4 


SEQUENTIAL DETECTION METHODS 


In Chapter 3 on fixed-sample-size methods, the underlying assumption was that in order 
to implement an FSS test, we must know the fixed-sample-size that fits the data. In ELINT 
applications, this is precisely the value that has the least a priori knowledge attached to it. We 
make assumptions regarding the bandwidth, center frequency, and modulations we desire to be 
captured, but the length of pulses is less well known. We only assume that the shortest pulse 
will be one sample long, and a longest pulse designation is only for the convenience of excluding 
CW carriers from our data collection. As such, there is no one “best” fixed-sample-size with 
which to accomplish detection. FSS methods perform poorly when the observation interval 
is wrong, by either reducing the effective SNR of received pulses, or performing suboptimally 
when the window length is too short. 

A sequential test attempts to rectify this situation. In this method, after each observation 
one of three decisions is made: accept hypothesis Hq, reject hypothesis Ho, or continue the test 
with another observation. As such, the length of the test n is a random variable, dependent 
on the data. For ELINT applications, this is particularly of interest as the testing procedure 
adapts itself to the data being received. Intuitively, we would assume that for a high SNR 
pulse, we may not need as many observations to reach a conclusion, but for low SNR pulses, 
provided that the pulse is of sufficient length, we could improve performance by allowing the 
test procedure more samples upon which to make its decision. 

The difference between FSS methods and sequential methods is seen most intuitively by 
noting the way in which each partitions the parameter space. In the last chapter, we denoted 
the critical region C as the area under which the acceptance of Hq was strongly preferred. In 
general, let us define an n-dimensional parameter space by R n . In FSS testing, according to 
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Figure 4.1 FSS partitioning of the N-dimensional parameter space into two mutually exclusive 
regions corresponding to the size and power of the test. 


the Neyman-Pearson theorem, for a fixed size a , we maximize the power of the test (3. Thus, 
we have split the parameter space into two mutually exclusive regions. For a fixed test of 
length N, R% is of size a and is the zone of preference for acceptance of H 0 , and R l N is the 
zone of preference for rejection of Hq. In general, it is a decision rule D which partitions the 
parameter space in this way by means of a sufficient statistic and critical values associated 
with the boundaries of each space. Graphically, we represent this process by a transformation 
D which takes the observation space Z and projects it into our divided parameter space to 
make a decision, as in Figure 4.1. Now, let us consider the partitioning of the parameter space 
under sequential procedures. Here, the space is divided into three mutually exclusive regions. 
Since the test length n is a random variable, the parameter space changes dimension with each 
successive sample, as do the shapes and sizes of the critical regions. We define another decision 
rule D* which takes the sample space and projects it into one of three regions. The first region 
is the zone of acceptance of Hq and is denoted by R^, which can be any size, as in the Neyman- 
Pearson formulation. Similarly, we have a zone of preference for the rejection of hypothesis 
Hq , R) n which, unlike FSS methods, can be made any size, since we now have another region 
which can cover all other possibilities. This region R m is the region of indifference, in which it 
is preferable neither to reject nor accept the hypothesis, but rather to continue observation, as 
we are yet unsure, based on the observations, what the correct course of action is. Graphically, 
we show the partitioning in Figure 4.2. The intuitive benefit of sequential testing is the ability 
to make the critical regions of any size, while leaving the remaining portion of the parameter 
space to a zone of indifference. In this paper, we will only consider tests which have probability 
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Figure 4.2 Sequential test partitioning of the parameter space into three mutually exclusive 
regions, the zone for preferable acceptance of Ho, the zone for preferable rejection of H 0 , and 
a zone of indifference where more observations are taken until absorption. 


1 of terminating. In other words, we may think of sequential methods as a way of adaptively 
computing the test size based on the data, which is a desirable property for ELINT. 


4.1 Sequential Probability Ratio Test 

So far, we have not said how these tests are constructed or what their properties are with 
respect to test length and SNR detection performance. In pulsed radar detection, the test 
length is important to the performance of our receiver. For example, if optimal sequential 
detection required 40 samples to detect on average, but on average we only expected to receive 
20 samples from each pulse, then we would conclude that despite its benefits, sequential testing 
does not perform to the required standard. Sequential testing methods are primarily the result 
of work by A. Wald [8] who is considered the originator of much of the statistical theory 
regarding sequential testing. However, the idea of sequential detection has been developed by 
several others including Page [9] and Lorden [10]. The basic theory behind sequential testing 
will be drawn from their work, and we will leave out some of the finer mathematical points 
not of particular interest to us, but will point out where more information regarding certain 
peculiarities can be found in these works. 

Let L n (9; Z\, Z 2 ,..., Z n ) be the likelihood of the observed samples Z n = [Z\, Z 2 ,..., Z n ] 
coming from a distribution with parameter 0. Then we could construct a likelihood ratio for 
testing the probability of a sequence of observed points to determine whether they were more 
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likely to come from one of two possible distributions, as in 


l n (0', Z n ) 


Ln{@'i Z r 


L n ( 0;Z n ) 

where the null hypothesis distribution parameter is assumed to be zero in the signal-absent 
case. This is not in general necessary, as signal absent may be a distribution whose parameter 
is not zero, but the results are general to any distribution with any parameter, provided that it 
is a complete family of distributions. The length of the sample set Z n is a random variable as 
in all sequential procedures. We now define our decision rule D* such that we continue taking 
observations whenever 


B <l n (9-Z n ) <A (4.2) 

for two positive constants A and B provided B < A. If 

ln{9 ; z„) < B (4.3) 

the test is terminated with the acceptance of Hq, and if 

ln{0; Z n ) > A (4.4) 

the test is terminated with the rejection of Hq. In this manner, we divide the parameter space 
into three mutually exclusive regions with regard to our preference for acceptance and rejection, 
as well as the zone of indifference. The likelihood ratio l n {0\ Z n ) can be expressed as the joint 
pdf of the samples of Z n as in 

“ ’ /(Zi,0)/(Z 2 ,0)---/(Z„,0) <4 ' 5) 

where every sample is i.i.d. with marginal pdf f(Zi,6). Suppose that we were to take the 
logarithm of the likelihood ratio to produce 

>o g Ue; z„) - log + log + • • • + log (4.6, 


We will denote the ith term of this sum as 


Zi = log 


f(Zi,0) 


1 0 ) 

We can reformulate the decision rule D* by replacing the constants A and B with log A and 
log B to achieve the same partitioning as above. If 

log B < zi + z 2 H- h z n < log A (4.8) 
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then continue the test with sample statistic z n+ 1 , but if 


zi+z 2 + -1- z n < log B (4.9) 

the null hypothesis is accepted, and if 

Zi + Z 2 -I -+ z n > log A (4.10) 

then the null hypothesis is rejected. This test is referred to as the sequential probability ratio 
test (SPRT) as developed by Wald [8]. 

4.2 Determining A and B 

As in constructing all statistical tests, we must determine the critical constants that ap¬ 
propriately size the parameter spaces f?°, R„, and R n . Consider first the rejection of the null 
hypothesis when 


l n {0,Zi n ) > A (4.11) 

Essentially we are saying that the joint probability measure f(Z\, 6 )f(Z 2 ,9), • • • , f{Z n ,9) is at 
least A times as large as the probability measure f(Zi, 0 )f(Z 2 , 0 ), • ■ ■ ,f{Z n , 0). We defined 
earlier the power of the test as /3 when H\ is true, which is the probability of rejecting the null 
hypothesis when it is false. However, there is associated with that the size of the test a which 
is the probability of rejecting the null hypothesis when it is true, or a false alarm. The joint 
probability measures above therefore state that for the hypothesis to be rejected, the power 
must be A times the false alarm probability. Thus, we generate the following inequality: 

A<- (4.12) 

a 

Similarly, for a case where the null hypothesis is accepted, 

l n (9,Z n )<B (4.13) 

The probability measure under the null pdf must be g times as large as that under the alter¬ 
native pdf. Associating the power and size of the test with this inequality we get 

B > £ (4.14) 

1 - a v ' 
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These inequalities yield the bounds of the constants, but we are interested more in determining 
their values for practical application.* For the above derivation we have found inequalities, but 
what is the outcome if we were to set the constants A and B to equalities in practice. It is 
sufficient to say in this instance that the determination of the exact power and size based on 
any two constants is quite tedious. Wald provides a method for this type of calculation, but also 
finds that for the strict equality, the size and power of the test do not change by an appreciable 
amount when the test is untruncated. The problem in determining the exact error bounds lies 
in what is called the “excess over boundaries” problem. In the above derivation, we assumed 
that the for the size and power of the test to have exact values a and /3, we achieved the bound 
A or B with strict equality. However, this might require a noninteger number of observations. 
When we are restricted to integer observations, as we always are, we may exceed the boundaries 
by a small amount, either decreasing the size or power of the test. Strict equality of the above 
constants to their derived values does yield one bound that is exact. That is, for a choice of 



(4.15) 


and 


B = 


1-P 

1 — a 


(4.16) 


either the size or power of the test will be exact, and the other will be decreased by a small 
amount due to the excess over the boundaries. 


4.3 Operating Characteristic Function 


In FSS methods, the probability of detection is a function of the size of the test and the 
limiting distribution of the likelihood statistic T. This is easy to find and simulate via math¬ 
ematical software because the length of the test is fixed. For sequential methods, however, in 
which the observation interval is itself a random variable, the distribution of the statistic is a 
function of the interval and the decision function. As such, asymptotic approximations to the 
actual performance of the test are the methods by which the tests are characterized. Suppose 

^The reader familiar with Wald’s work may notice the difference between Wald’s bounds of A = and 

jj Oc 

& = fr° m those presented here. There is no difference, however, except in notation. Earlier, we presented 
the power of the test as f3, but in Wald’s derivation the power of the test is 1 - (3 with being the probability 
of a Type II error. Thus, we have used /3 here to denote power to be consistent with earlier sections, but the 
bounds remain the same regardless of the different formulation. 
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we have a function C(6) which is the conditional probability of accepting hypothesis H 0 at the 
parameter point 9. Then for no signal, or 0 = 0, C(0) = 1 - a and is otherwise a function of 
the parameter 0 , C(0) = 1 - 0(0). Thus, the probability of detection for 0, or the power of the 
test, is 


m = i-m 


(4.17) 


We may specify one value of (3(0) at a parameter point 0 = 9 d to yield the bounds A and B. 
To determine C(9), consider the statistic 


£(fli;Z n ) l* w 

L(o-,z n )_ 


(4.18) 


where, for every value 0 , the value of h(0) is determined such that the expected value of 
Equation (4.18) is 1, as in 


L(9 i;Z n ) 1 ft(tf) 

L( 0;Z n ). 


L(0] Z n ) dZ n — 1 


(4.19) 


which was derived by Wald [8]. It follows that the integrand of Equation (4.19) is a distribution 
of Z n which we denote by 


1 (Zn) “ mzT) J i(S ' Zn) 


Let H denote the hypothesis that 


m Zn) = f(Z l ,9)f(Z 2 ,0)--- f(Z n , 9) = f( Z n , 0 ) 


(4.20) 


(4.21) 


is the true distribution of Z n , and H* the hypothesis that /*(Z n ) of Equation (4.20) is the 
distribution of Z n . Consider a SPRT which continues taking observations when 


Bm < < A m 

J v) 


(4.22) 


accepts H* when the ratio is greater than or equal to A h ^ and accepts H when the ratio is 
less than or equal to B h ^ e \ From Equation (4.20), 


/*(Z n ) r m;Z n ) l fe W 
/(Z n , 9) [t(0; Z„) 


so that Equation (4.22) can be rewritten as 


L^) 

L(0;Z n ) ^ 


(4.23) 


(4.24) 
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which is the standard SPRT derived in Section 4.1. 

If the test between H* and H results in the acceptance of H*, then Equation (4.24) implies 
the acceptance of Hi. Likewise, the acceptance of H implies Hq. It follows that £(0), the 
probability of accepting H 0 , given 6 = 0, is the same as £(0), the probability of accepting H 
when /*(Z n ) is the true distribution. To calculate £(0), let a' and /?' be the size and power of 
the test of H* versus H, respectively. It follows that 

A h ^ = ^ (4.25) 

and 

Bh{9) = (4 - 26) 

Solving for a' and noting that 


C{6) = 1 -a' (4.27) 

the operating characteristic function (OCF) is given by 

A h ( e ) _ i 

= AH6) - bW) ( 4 - 28 ) 

when we neglect excess over boundaries as before. 


4.4 Average Sample Number 

Given the OCF derived in Section 4.3, we can use that information to determine the average 
number of samples required for detection depending on a particular implementation of the 
SPRT. Let 


Z n = log 1(0; Z n ) = zi + z 2 H- + z n (4.29) 

as in the log-likelihood formulation of the SPRT in Section 4.1. The test procedure is as usual 
for a SPRT: reject H 0 when Z n > log A, accept H 0 when Z n < log B, and continue the test 
when log B < Z n < log A. If we calculate the average value of Z n at a test that has terminated 
at length n, and we neglect excess over boundaries, then the expectation is given approximately 
by log B times the probability of accepting Ho plus log A times the probability of rejecting H 0 , 

E{Z n , 9) = C(9) log B + [1 - C{9)} log A (4.30) 
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(4.31) 


From the equation for Z n , 

E{Z n ) = E(z\ + Z 2 + ■ • • + z n ) = E(n)E(zi ) 

for independent z { , from i = 1,... ,n, giving E(z) = E(z x ) = E(z 2 ) = ■■■ = E{z n ). Thus, we 
denote the average sample number (ASN), n — E(n) as 

_ , C(6) log R + [1 — C(Q)\ log A ' x 

n = E(n) = & - ~E{z) (4.32) 

Given the OCF of a particular sequential test, we can use that to find the average number 
of samples that will be required for a particular decision, either the acceptance or rejection 
of Ho- This will be useful to us later in the design of sequential tests to achieve the desired 
performance, and as a way of characterizing their effectiveness versus other sequential tests or 
FSS tests. 

4.5 Cumulative Sum Tests 

Suppose that we are only interested in a one-sided hypothesis test; for example, we only 
wish to detect the acceptance of H x regardless of the times when Ho is accepted. Or, consider 
the possibility that we want a test for multiple hypotheses. In these cases, one implementation 
of the SPRT becomes useful, the cumulative sum test, or CUSUM test. The CUSUM was 
originally discovered during manufacturing processes when the quality of the output slowly 
drifted with time, and FSS estimators were unable to detect the slow change. The CUSUM 
relies on the properties of the logarithm in the log form of the SPRT in that for a single point 
that is more likely noise than signal, the likelihood ratio will be less than one, resulting in a 
negative log. For the contrary case, where the point is more likely to come from signal than 
noise, the likelihood ratio is positive. Thus, for periods in which noise is the predominating 
factor, the slope of a cumulative sum of each likelihood ratio point 

Z n — z\ + z 2 + • • • + z n (4.33) 

will be negative, and a positive slope will result from a signal present. Thus, let us only concern 
ourselves with the bound associated with rejection of Ho , namely when 

Z n > log A (4.34) 
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We denote the cumulative sum by the notation CUSUM(k,i) where 

k 

CUSUM(i , k)='£z t = Z k -Z i (4.35) 

t—i 

In the case of the cumulative sum, there is never an acceptance of Hq; the test merely continues 
until there is some rejection of H 0 . To determine the stopping time of this statistic, we cite 
the work of Lorden [10]. Lorden determined that any difference between the current CUSUM 
value and any previous value that exceeds Wald’s threshold log A would result in the same 
performance as that of the standard SPRT. Thus, the stopping point is the first point such that 

(4.36) 

Using this method, we achieve the same results as the standard SPRT, but in some cases it 
results in a simplified implementation, when the acceptance of H 0 is not important. As an 
example, in Figure 4.3(a),(b) we show the CUSUM as it would be formulated for the detection 
of a shift in the mean of a Gaussian process. We insert a shift in mean from the 1000th to the 
2000th sample. During that interval, the CUSUM generates a positive slope. In Figure 4.3(c), 
we show graphically the interpretation of the stopping rule of Lorden, Equation (4.36). 

One other advantage of the CUSUM rule is that we can use it as a multihypothesis testing 
statistic. Consider the problem of estimating the time of arrival (TOA) of a pulse as well as 
the time of departure (TOD) or end of the pulse. In the standard SPRT, there is no good way 
of knowing the time of return to the previous statistics. However, for the CUSUM test, we can 
find the bound easily, by looking for a negative slope that satisfies Wald’s bounds. Thus, as in 
Figure 4.3(c), we find the first point such that 


n i = inf 


k > 1 : max CUSUM(i,k) > log A 

Ki<k 


n 2 = inf 


t > N\ 


max CUSUM(i, k) < — log A 

Ni <i<k 


(4.37) 


The cumulative sum test is important for ELINT applications as it gives us both the beginning 
and end of the pulse with similar error bounds. It also alleviates the problem of the SPRT 
alarming multiple times within a pulse. Since the SPRT resets itself at the acceptance or 
rejection of H 0 , in the duration of a long pulse, we may expect the SPRT to alarm multiple 
times, which causes needless extra complexity to other processes which rely on accurate TOA 
times. 
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(a) Gaussian Shift in Mean Process 


(b) CUSUM 



(c) Stopping Rules 


Figure 4.3 (a),(b) Example of CUSUM output for a Gaussian shift in mean process from the 
1000th to the 2000th sample, (c) Graphical interpretation of the stopping rules. 


4.6 Optimal SPRT Tests for Detection without Interference 

In this and all future sections, we will address the SPRT in general, but the results are 
applicable as well to the cumulative sum implementation of the SPRT, and the performance 
characteristics in terms of detection probability and average sample size are identical. We now 
address the problem of finding a SPRT test suitable to a particular set of statistics. This 
involves mainly the derivation of the appropriate log-likelihood ratio statistic, and setting the 
power of the test at a certain parameter point Oj, so that the test achieves desirable results on 
average. 
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4.6.1 Gaussian SPRT 


As an example that will be influential later in this paper, let us derive the log-likelihood 
ratio statistic in the case of a Gaussian process where we wish to detect a shift in the mean, 
assuming the variance of the noise to be constant. The likelihood ratio test for a Gaussian 
process is of the form 


i jZi-hr 


__ e 2 V 


(4.38) 


where the log-likelihood is given by 


Zi = ]Qgl(Z i] e 1 ,Oo) = £(z i -9o-£ 


(4.39) 


where <j> = 6 i — 9 0 and after a series of algebraic manipulations. The test is defined as the sum 
of these variables, Z n = z\+ z^-\ -+ z n , as in 

£(*-*-?) (4-40) 


= (4.40) 

This test is particularly attractive to FPGA implementation due to its ease of implementa¬ 
tion, requiring only an estimate of the variance and the mean, as well as additions and two 
multiplications. This results in relatively little complexity for implementation. 


4.6.2 Chi-Squared SPRT 

From the section on the polyphase filter bank, we know that our actual statistics are chi- 
squared random variables. In order to implement the optimal sequential test, we must first 
derive the log-likelihood ratio as in the case of the Gaussian test. Thus, for a parameter 
which in actual cases is the noncentrality parameter A^, we find that the likelihood ratio is of 


the form 




(4.41) 


For this derivation we adopt the noncentral chi-square distribution formula developed by Fisher 

[ 11 ] 


P(*>*) = |e x t X I 0 (VXx) 


(4.42) 
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where Iq(x) is the modified Bessel function of the zeroth order. Taking the logarithm and 
canceling terms, we get 

Zi = ~~2 + l0g (4.43) 

This SPRT statistic is optimal in the sense of fitting the expected statistics of the output of the 
polyphase filter bank. However, for FPGA implementation, the results are nearly impossible to 
implement because of both the logarithm and the modified Bessel function. Both of these would 
have to be stored in lookup tables that are expensive in resources to implement. However, let us 
look at reasonable approximations to the optimal solution. Suppose we expand the logarithmic 
Bessel function in a power series contingent upon a weak signal in which A << 1, which gives 
us 

* = -y + T Zi-^ + oWzf) (4.44) 

It may be noticed that the third term in Equation (4.44) contributes to the bias of Zi, and 
without that contribution the average sample number approaches infinity for hypothesis Hq. 
Keeping the terms up to second order results in the need to square each output value. This is also 
unacceptable in FPGA implementation because of the excessive number of bits that are needed 
for the computation. A conservative estimate on the output of the filter bank would be 16 bits 
for an accurate representation of the dynamic range of the input to the ADC. Implementing 
this form of the SPRT would require us to nominally keep 33 bits for the SPRT, which is 
far too many to provide efficient implementation. The assumption that A << 1 is actually 
a very reasonable approximation for FPGAs, as a common method is to use fractional two’s- 
complement representation internally. In this method, the numbers are stored in a fractional 
context with the maximum absolute value of 1. Thus, any noncentrality parameter that does 
not overflow the ADC would necessarily be less than 1, and much less than 1 in the case of all 
but the highest SNR pulses where the approximation need not be so accurate as those pulses 
are easy to detect. 

Let us continue approximations that may yield a useful result. From the properties of the 
chi-squared distribution it can be shown that 

E(Zi) = 2 + A (4.45) 
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and 


E(Zf) = 8 + 8A + A 2 


(4.46) 


Thus, for a received signal with noncentrality parameter A, the expected value of the test 
statistic is 

E( Zi \\) * ^ ^ (4.47) 

when terms higher than second order are ignored due to the negligible effects when X d « 1 
and A << 1. The second term in the approximate expectation is the average bias contribution 
of the squared Z x term in the small signal expansion. Suppose now that we keep only the first 
two terms of the small signal power series and the average contribution of the bias of the third 
term. It may then be shown that the test will terminate with probability 1; however, average 
sample numbers with the approximate test may be longer than the optimal implementation 
using the true bias term or even an exact modified Bessel function. We now have that 


Ad _ A^ A dZj 

2 8 4 


(4.48) 


For FPGAs, this approximation is attractive because we have made no assumptions that will in 
general be violated for our standard operating conditions, and we have reduced the complexity 
of the implementation to a manageable level. But the question is, what does this approximation 
tell us? Let us denote the mean of a chi-squared random variable with two degrees of freedom, 
as at the output of the polyphase filter bank in the case of noise, as fi x = 2 and the variance of 
the same variable as cr£ = 4. Rewriting Equation (4.48) in an equivalent form we have 



which is exactly equivalent to the Gaussian SPRT statistic of Equation (4.39). Thus, through 
the small signal approximations we have concluded that the limiting distribution of this set of 
statistics is Gaussian. This is an important conclusion, as we have derived an optimal small 
signal chi-squared detector that is equivalent in performance to a Gaussian detector as would 
be implemented in a coherent detector where the phase of the signal is known. As such, we 
expect that the test will achieve results that are comparable to matched filter outputs in the 
case of small signals and fractional two’s-complement binary arithmetic. 
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To characterize the performance of this test, we now derive the OCF of the implementation 
of the SPRT in Equation (4.49) with the methods described in Section 4.3. As such, we wish 
to find h( A) such that 

<4 - 5o) 

Using the small signal approximations above, it can be shown that the solution to this equation 
is approximately 


h{ A)«l —2^ (4.51) 

which, it is interesting to note, is the same h{9) that can be derived from the Gaussian SPRT 
[12], Thus the OCF for this test is given by 


£(A) = 


_ i 

.41—2A/A d _ gl-2A/A d 


(4.52) 


With the OCF, we can also derive the average sample number using the expected value of zt 
from Equation (4.47). As in Equation (4.32), 


= _ £(A) log B+[l- £(A)] log A 

4 8 


(4.53) 


4.7 Performance and Design of SPRT Tests 


For general SPRT tests, we have so far shown that the performance is based both on the 
critical bounds A and B as well as the design parameter 6 d . In the case of the chi-squared 
SPRT, the design SNR, represented by the design noncentrality parameter X d , is an important 
consideration, as well as assigning the power of the test at that SNR. Usually, the probabil¬ 
ity of false alarm is fixed in any application by practical considerations of data volume. In 
this project, the design value was a FAR of 10 -3 and we will use that value for subsequent 
calculations, although the design methods apply to any FAR. It can be seen from the section 
on determining the critical constants that one of the advantages of the SPRT is the ability to 
design an arbitrarily powerful test at any design point. We can pick any (3 at any design point 
A d without restriction except that 0 < [3 < 1. 

However, there are several considerations that the designer should bear in mind when con¬ 
structing these tests. High power intuitively means that we are accepting a high average sample 
number as a trade-off. In fact, it will be seen that the ASN of a SPRT is highly nonlinear, and 



in some respects can only be a component of designing tests by trial and error. Let us outline a 
general method of design and provide some examples that will give insight into the important 
considerations when designing these tests. 

In practical applications, it is often important that we achieve detection in as little time as 
possible, and it has been shown by Wald that the SPRT minimizes the time to detection of any 
test, when the optimal likelihood ratio is used [8]. Often, we must sacrifice test power to achieve 
a reasonable detection time. In the case of FSS methods, it was shown that overestimating the 
observation interval resulted in reduced SNR. In the SPRT, the effect is similar, except that 
the effect is not on SNR, but on time to detection. If the time to detection is on the average 
too long, then we run the risk of coming to the end of a pulse before detection. This problem 
motivates the first step in SPRT test design. 

Let us define a general SPRT test S that is implemented in the standard method. For 
each /3 on the range from 0 < 3 < 1 we can calculate the average sample number associated 
with a received signal at the design parameter A^. For increasing power, we expect increasing 
average sample numbers. Let us define an average sample number n which is the desired time 
to detection for a pulse at SNR commensurate with the design point. From the curve of power 
versus average sample number, we can read off the maximum power that can be achieved at 
the design ASN of n^. To verify that system performance is adequate, we compute the OCF 
and ASN for the SPRT test S given the power found from the design ASN, and verify that the 
system responds appropriately at all other SNRs. 

Consider the example of the chi-squared SPRT with a design parameter of -5-dB and a 
design ASN of 8. In other words, we would like to design a test that alarms on the output of 
the polyphase filter bank and achieves on average a detection time of 8 samples for a -5-dB 
wideband SNR pulse. This is a reasonable design goal, as it amounts to a detection of a -5-dB 
pulse that is 2.5 /us or greater. We compute the ASN as a function of test power and note, 
from Figure 4.4(a), that we can generate a test power of .61 for a design ASN of 8. If we 
then compute the OCF and ASN of this test based upon a power of .61 and a FAR of 10~ 3 
at a design SNR of-5-dB, we can verify that the performance of the test at all other SNRs is 
reasonable. From Figure 4.4(c), we see that the average sample number peaks at slightly less 
than -5-dB at approximately 9 samples, which is not unreasonable at lower SNR. Thus we may 
conclude that this design meets the specifications we imposed at the outset. 
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Now, let us consider another design process where we set the design SNR slightly lower, at 
-10-dB and a design ASN of 20. If we again compute the ASN as a function of test power, we 
see in Figure 4.5(a) that we can achieve a maximum power of approximately .26 for these design 
constraints. If we now plot the probability of detection derived from the OCF and the ASN 
associated with the test, we can see in Figure 4.5(c) that the design provides an ASN of less 
than 20 at all other SNRs, which seems reasonable. However, let us compare the performance of 
this design at -5-dB with the design previously where the design parameter was -5-dB. Looking 
at the probability of detection for the second design in Figure 4.5(b), we see that at -5-dB we 
achieve only a power of .78 at the expense of an ASN that is nearly twice that of the first 
design. If we set out with an ASN of 15 in the first design, we could achieve greater than a 95% 
probability of detection. Thus, for a lower SNR design parameter we have sacrificed heavily 
the optimal performance that could have been achieved at higher SNRs. 

This example illustrates an important point. Although we may achieve the design goals of 
a particular design ASN and design SNR, we may be unduly sacrificing performance at other 
SNRs if we choose incorrectly the design parameters. It is important to design SPRT tests of 
this sort through trial and error to achieve acceptable performance at all SNRs. In general there 
is no optimal way to design these tests so that all SNRs are arbitrarily powerful with acceptable 
ASNs. Therefore, use of these design methods with a range of design parameters is useful to 
pick a design parameter that is most favorable in performance over the entire detection range. 


4.8 Optimal SPRT Tests for Detection with Interference 


As in Chapter 3, we consider the detection of pulses in the presence of CW interference. We 
saw in the last chapter that in the presence of interference the noise pdf assumes a noncentral 
chi-square distribution, giving the likelihood ratio as 


\d) 


<r 2 xj( M 

<7 2 X2 2 ( a o) 


(4.54) 


where A rf is the design parameter, and Ao is the noncentrality parameter of the interference. If 
we use Fisher’s formula for the y' 2 2 (A) pdf as we did in the first derivation, it can be shown that 
the optimal log-likelihood SPRT test is given by 


Zi = 


~ A o , (hW^d,Zi) ^ 


(4.55) 
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where Iq(x ) is again the modified Bessel function of order 0. Let us apply the same small-signal 
approximations to Iq(x) as before, keeping the first two terms of the power series expansion 
and the average bias contribution of the third term. We can then show that the small-signal 
approximation yields 

*=^( Zi -(A„ + 2)-V^) (4.56) 

which is exactly the Gaussian SPRT with the mean of the noise correctly given as = 2 + Ao 
and the variance as cr£, = 4. Thus, we may conclude that the Gaussian approximation is 
equally valid in the case of CW interference and we need not change critical bounds in this case 
to achieve similar performance. 

4.9 Comparison of Performance of SPRT versus FSS Methods 

Suppose that we want to characterize the performance of the SPRT methods presented in 
this chapter with the FSS methods presented in the last. Let us first consider the average 
sample size to detection, and see what the trade-offs between each implementation are. We 
begin by assuming that we are using a 32-sample window for Zn x and a 64-sample observation 
interval for Zn 0 . We compute the probability of detection using the F-test for detection of 
signals with out CW interference. If we then compute the ASN for an SPRT that is designed 
to achieve the same detection probability curve, we can compare the ASN of the SPRT to 
detection versus the 32 points required for the FSS test. In Figure 4.6 we see that, while the 
FSS methods require 32 samples to achieve the same detection probability, the SPRT requires 
only a maximum of 7 and in some cases as little as 1 sample. 

Now, for a similar sample number, let us compare detection probability. Consider the two 
examples presented in Section 4.7, in which we designed an SPRT based upon an ASN design 
value and a design SNR. If we design an SPRT in this fashion, and compare its detection 
probability with an FSS test with sample size equal to the maximum ASN of the SPRT, what 
is the difference in detection probability? In Figure 4.7 we see that in most cases the SPRT 
outperforms the FSS tests significantly, except for a small region in the SPRT designed for 
-10-dB. Thus, we conclude that in all cases, we can design an SPRT that will be superior to 
the optimal FSS test in both ASN and probability of detection. 
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4.10 FPGA Implementation 


Under the small-signal assumptions, the chi-square SPRT is easily implemented on FPGAs. 
If we replace n x > in the CW interference derivation of the optimal SPRT, Equation 4.56, with 
its ML estimate, the sample mean, we can implement one detector for both interference and 
noninterference. Likewise, we can replace the variance with an estimate of the variance from the 
data to account for other inconsistencies in the true distribution of the observed data. In this 
way, we construct a robust detector that is capable of achieving desired performance in a range 
of adverse conditions, including interference, nonstationarity, and fixed-precision effects. In 
addition, this detector fits the statistics well with no unreasonable approximations and results 
in a low complexity implementation. Thus, we conclude that for ELINT applications using 
subband channelization based on the DFT, the SPRT is a superior choice for detection over 
FSS methods. 
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(a) ASN vs. Test Power 



(b) Probability of Detection ( c ) ASN 


Figure 4.4 (a) ASN of the SPRT test as a function of test power at a design parameter of 
-5-dB SNR. The dotted lines show the maximum test power that can be achieved at a design 
ASN of 8. (b) Probability of detection and (c) ASN for a design parameter of -5 dB and a 
design ASN of 8. 
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(a) ASN vs. Test Power 
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(b) Probability of Detection (c) ASN 

Figure 4.5 (a) ASN of the SPRT test as a function of test power at a design parameter of 
-10-dB SNR. The dotted lines show the maximum test power that can be achieved at a design 
ASN of 20. (b) Probability of detection and (c) ASN for a design parameter of -10-dB and a 
design ASN of 20. 
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(a) Probability of Detection for both FSS (b) ASN and FSS Sample Size 

and SPRT Tests 

Figure 4.6 (a) Probability of detection for an FSS test of observations No = 64 and Ni = 32. 
The SPRT is designed to achieve the same detection probability, (b) ASN of the SPRT and 
sample size of the FSS test to achieve the same detection probability. 


Companion of Probability of Detaction for FSS vs SPRT Companion of Probability of Detection for FSS Vi SPRT 



Wideband SNR Wideband SNR 


(a) -5-dB Design SNR and 8 Sample Design (b) -10-dB Design SNR and 20 Sample De- 

ASN sign ASN 

Figure 4.7 (a) Probability of detection for SPRT and FSS methods based on a design SNR 
of-5-dB and design ASN of 8 for the SPRT, and a FSS sample size of 9. (b)(a) Probability of 
detection for SPRT and FSS methods based on a design SNR of -10-dB and design ASN of 20 
for the SPRT, and a FSS sample size of 21. 
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CHAPTER 5 


CENTER FREQUENCY AND BANDWIDTH ESTIMATION 


As a further requirement of this project, we are charged with determining the bandwidth 
and center frequency of a detected pulse for output to further processing elements that will 
digitally modulate the signal to baseband and adaptively decimate the signal to the matche 1 
bandwidth estimate of the modulation. If possible, we would like to exploit our estimators used 
for detection to implement bandwidth calculations. 

There are important design goals and limitations to the bandwidth measurement scheme. 
We desire to provide the system with dehop/dechirp functionality. In other words, we would 
like to follow the center frequency of the pulse in the case that it is linearly modulated as in 
the case of LFMOP, or in the case that it takes randomly assigned frequencies as in the case of 
a hopped NOMOP signal. The limitation imposed upon this functionality is that we can only 
update the digital baseband modulator on a 1 (is basis. 

5.1 FSS Detection Bandwidth 

Given a 1 /us basis and a decimation rate of 32 in this project after the DFT, the best 
we can do at an integer sample size is 4 subband samples, representing a 1.28 fxs update rate 
in the undecimated data stream. Thus, our goal is to update both the bandwidth and center 
frequency every four subband points after detection, and make that information available for 
further processing. For the FSS detection methods of Chapter 3, the bandwidth estimation is 
a further application of FSS methods. Let us define our TOA detection point as T. We would 
like to retain only statistically significant energy in each channel at a specified confidence level 
OiBW- For this, it seems reasonable to implement the F-test based on four observations from the 
pulse against Nq observations of the noise, just as in the FSS scheme presented in Chapter 3. 
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Thus we calculate the statistic 


1 \~^K(x)+3 y 
4 2^i=K(x) 

7 % 


(5.1) 


where K(x) = [T, T + 4, T + 8,... ] encompassing every fourth point from the TOA to the end 
of the pulse, and Y^n 0 indicates the sum of all points in the window N 0 for the calculation 
of the noise background. Then, if the statistic 


— F8 t 2No;aBW (5.2) 

we may conclude that there is statistically significant energy in that channel over 1.28 (is 
to be included in the significant bandwidth. Similarly, in the case that detection is made in 
interference, the critical bound is given by the noncentral F test. 


5.2 SPRT Detection Bandwidth 

In the case that the detector is based upon SPRT methods, as is preferable, we would like 
to exploit that statistic to perform bandwidth estimation without having to recompute FSS 
statistics. Under the small-signal model, we have an implementation which takes the current 
point and subtracts an estimate of the mean along with an additional factor based on the 
design SNR. If we implement the SPRT in the cumulative sum fashion, we notice that this 
particular implementation is nothing more than an approximate average short-time Fourier 
transform (STFT) as a result of the polyphase filter bank. The averaged STFT is referred to as 
a periodogram and is a common method for reducing the variance of the DFT of a stationary 
process. Thus, to determine significant bandwidth, we need only see which averaged CUSUM 
point exceeds the background energy level. Consider the statistic 

U = CUSUM(K(x ) + 3,1) - CUSUM(K{x ), 1) (5.3) 

where K(x) is the same as for FSS estimation. Over a range of four samples, if the CUSUM 
detects a more probable signal than noise, the slope will be positive on average, and if it is more 
likely noise than signal, it will be negative. Thus, a reasonable bandwidth estimation scheme 
would be to construct these statistics with a critical bound of 0. If 

U > 0 (5.4) 
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then we conclude that there is significant energy in that channel. In this case, since we are using 
a mean estimate for the noncentrality parameter of the interference in the case of detection in 
interference, it may be necessary to subtract an additional confidence factor from the statistic 
U to account for the variance in the estimate of the background noncentrality parameter. 

5.3 Advanced Bandwidth Estimation Methods 

It is possible to apply advanced methods to the bandwidth estimation problem to determine 
modulation and other interesting factors about its behavior, but these methods are too complex 
to implement in today’s FPGAs. We provide a note here about some interesting methods that 
may be possible in the near future when this amount of complexity can be implemented. Under 
the small-signal model, we can roughly approximate the statistics as being Gaussian. T his 
approximation allows a whole body of theory dealing with the analysis of variance (ANOVA) 
to be exploited for bandwidth detection. Under these models, we would be able to say, for 
example, whether there is significant energy in any number of channels simultaneously, or 
whether interaction between these channels at any stage of detection is significant. This is a 
powerful method for determining modulation. We can use methods by Scheffe [5] to construct 
confidence intervals for energy contained between and within channels, as well as models which 
would explain the interaction over time between channel background noise and the signal energy. 

5.4 Excluding Pulse on Pulse Detection 

In general, there are two ways of determining the bandwidth of a signal: one would be 
to take the total bandwidth encompassing any channel with significant energy, and the other 
would be to take only the bandwidth that is coherently distributed, i.e., occupies adjacent 
channels. We define a pulse-on-pulse situation as that in which there are two pulses occurring 
simultaneously with at least one band of separation between their effective bandwidths. In 
some cases, we would prefer to detect and estimate the bandwidth of only the highest-SNR 
pulse, eliminating the other through the adaptive filtering process. 

To eliminate a pulse-on-pulse bandwidth expansion from occurring, we would like to detect 
only coherent bandwidth, where the energy is distributed among adjacent bands. Therefore, 
a logical approach would be to start with the channel of highest energy and look out in both 
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(a) Normal Pulse Bandwidth 



Start from Maximum of STFT 
Periodogram and look out for Bandwidth. 
Stop at first negative STFT Band. 


“i - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 - 1 “ 

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 

Bands 


(b) Pulse-on-Pulse Bandwidth 


Figure 5.1 Bandwidth calculations excluding pulse-on-pulse occurrences. 


directions, towards Nyquist and dc, until the statistic U has a value less than zero. Thus, in the 
case that pulses coincide and are separated in bandwidth by at least one band, we can alarm 
on the highest-SNR signal. We can see in Figure 5.1(b) that the simultaneous pulse is ignored 
and the bandwidth of the highest-SNR pulse is retained. 

5.5 Pulse-on-Pulse Bandwidth Detection Problems 

While the pulse-on-pulse bandwidth calculations exclude typical simultaneous pulse matched 
bandwidth problems, they also cause problems with the bandwidth detection of certain kinds 
of MOP pulses in which the energy is not coherently distributed. As an example, we show a 
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STFT Periodogram for PSKMOP Pulse 



Figure 5.2 STFT from a PSKMOP pulse. 


typical periodogram from a PSKMOP pulse. In Figure 5.2 we see that although no pulse-on- 
pulse behavior is occurring, the modulation bandwidth will not be correctly identified given 
the noncoherent distribution of energy among the subbands. Thus, we conclude that a better 
strategy for capturing the true modulation bandwidth is to choose a bandwidth that covers 
all of the significant energy in any subband, which we call the min/max bandwidth estima¬ 
tion scheme. In Figure 5.3 we use the min/max bandwidth estimation to capture the true 
bandwidth of the modulation. In the case of a true pulse-on-pulse situation, the result would 
be the detecting and encompassing of both pulses in the matched bandwidth. Either of the 
methods discussed for measuring the bandwidth may cause problems depending on the nature 
of the post-detection processing. Thus, only an application specific choice can be made be¬ 
tween the merits of excluding pulse-on-pulse bandwidth at the expense of degraded PSKMOP 
performance. 

5.6 Decimation Rate Stabilization 

Because of FPGA implementation and our high input clock rate, the ability to filter the 
signal to matched bandwidth and then decimate causes important gains in FPGA resource usage 
as well as the volume of data produced. An adaptive decimator can achieve good performance 
at matching the pulse bandwidth; however, it is clear that switching the decimation rate in 
the middle of the pulse will cause unwanted transients and a differing time-frequency scaling 
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Bandwidth 



Bands 

Figure 5.3 Bandwidth calculation using the min/max bandwidth. 

throughout the pulse. As such, though we are capable of achieving a 1.28 n s update rate on the 
tuning of the adaptive decimator, it is unwise to change decimation rates inside a pulse. Thus, 
we seek some method to calculate a “best” bandwidth that encompasses the maximum amount 
of energy in the pulse. To do this effectively, the best method would be to take snapshots of 
the bandwidth of the pulse every 1.28 fis and choose the greatest bandwidth to be sure to 
capture all the energy over the length of the pulse. Because of FPGA memory limitations, we 
cannot, in the event of long pulse widths, buffer the entire pulse data. Thus, beginning with 
the TOA of the pulse AT(1), we must calculate a new bandwidth every fourth point using the 
methods of pulse-on-pulse estimation or min/max estimation. Since there will be some variance 
in the estimate, we would like to take multiple observations of the bandwidth and choose the 
largest. We define X B w as the number of bandwidth measurements to be taken before choosing 
a bandwidth tune word for the set-on receiver. In the event of long pulse detection, as X B w 
grows, the stability of the bandwidth measurement increases. 

5.7 Center Frequency Estimation 

After we have calculated our stable bandwidth over the duration of the pulse, what re ma ins 
for a full tune word to be output to the set-on receiver is the center frequency. 

A desirable property of any set-on receiver/detector pair would be the ability to dechirp an 
LFMOP signal and dehop a NOMOP signal. In the case of this algorithm, we implement this 
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Frequency 


Figure 5.4 Non-baseband-modulated asymmetric bandwidth. 

functionality by updating the center frequency on a 1.28 /zs basis. Once we have the stable 
bandwidth from Xbw observations of the bandwidth, we return to the TOA point generated by 
the SPRT and apply that stable bandwidth to finding the center frequency every fourth point. 
By allowing the center frequency to be updated on a 1.28 /zs basis, we provide for the ability 
to dechirp long-duration LFMOP signals. For short-duration LFMOP, less than the update 
rate, the full bandwidth is captured, with no dechirping effect. For dehop capability, no center 
frequency update is needed over a NOMOP pulse as the stable bandwidth should cause the 
center frequency to remain constant throughout a single pulse. In a pulse train, the separate 
detection of each pulse will provide a new center frequency in each case, resulting in dehopped 
output. 


5.8 Center Frequency Location 

In general, there are two methods for placing the center frequency given a stable bandwidth. 
One would be to place the center frequency at the highest energy level subband. The second 
would be to place the center frequency at the middle of the bandwidth that covers the sig nifi cant 
energy. To clarify, we offer this example, where there are the two options above. In Figure 5.4, 
we show a pulse situation with an asymmetric bandwidth, with the two options for center 
frequency calculation. Here, we assume the ability to adjust the bandwidth to cover the coherent 
energy spectrum. Choosing center frequency and bandwidth #1, as in Figure 5.5, we see that 
we introduce noise into the matched bandwidth filtering because of the asymmetric bandwidth. 
If we choose center frequency #2, as in Figure 5.6, we realize a true matched bandwidth process 




, —-I-- 

-7t/M Jt/M 

Figure 5.5 Excess noise introduced to matched bandwidth for center frequency #1. 


Center 

Frequency #2 



—7l/M k/M 


Figure 5.6 True matched bandwidth for center frequency #2. 

by optimally filtering the significant bandwidth. The only error introduced in this method would 
be in the case where the true center frequency is placed on the filter transition band edge, which 
would result in some attenuation of that center frequency due to nonideal low-pass filtering. 

5.9 Center Frequency Calculation 

Given the above choices for center frequency, we choose to locate our center frequency at 
the middle of the stable bandwidth for optimal matched bandwidth filtering. Once we obtain 
the stable bandwidth estimate, we must return to the TOA location if(l) and begin a center 
frequency calculation every fourth point. We use a method similar to bandwidth estimation by 
first calculating the statistic U. After we generate the periodogram, we would like to maximize 
the energy covered by the stable bandwidth. Given the stable bandwidth and the 27 valid 
steps of our baseband modulator, we realize some number of allowable center frequencies for 
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Subband 

Figure 5.7 Allowable tuning steps for center frequency 


the bandwidth to be contained between the dc and Nyquist frequency. By summing the energy 
contained in the periodogram under the stable bandwidth at each step, we simply choose the 
center frequency that generates the maximum energy. This results in a center frequency that 
optimally covers the most significant energy at each 1.28 /j,s update. 

5.10 FPGA Implementation 

Using the statistics already derived in either the FSS or SPRT methods provides sim¬ 
ple methods of bandwidth and center frequency estimation. Each method involves very little 
computational complexity. The advanced ANOVA techniques would require significantly more 
computational complexity, possibly involving a further transformation of the data, and squaring 
operations. These methods are unsuitable for implementation in this project at this point; how¬ 
ever, they may be useful in future applications on more advanced FPGAs or in a floating-point 
format. 
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